Here we replicate the paper’s work by bootstrapping the spectral density of the autoregressive process \[X_t = 0.5X_{t-1} - 0.6X_{t-2} + 0.3X_{t-3} - 0.4X_{t-4} + 0.2X_{t-5} + N_t\] using the given bootstrap algorithm:
set.seed(505)
x1 = arima.sim(n=256, list(ar=c(0.5, -0.6, 0.3, -0.4, 0.2)))
tsplot(x1, col=4)
mvspec(x1, col=4)
mvspec(x1, spans=c(2,2), col=4)
mvspec(x1, spans=c(2,2), log="yes", col=2)
spec = arma.spec(ar=c(0.5, -0.6, 0.3, -0.4, 0.2), main="Autoregression")
Usually, we can only create a confidence interval for our spectral
density by using mvspec with the log option.
The confidence interval is based on a Chi-square distribution with only
2 degress of freedom. Thus, the interval is too wide and of little use.
Therefore, in this work, we will try using a bootstrap approach to
create more accurate spectral density estimates.
The following step is to extract the size and frequency of the data.
xfreq = frequency(x1)
n = length(x1)
nspec = floor(n/2)
freq = seq(from = xfreq/n, by = xfreq/n, length = nspec) # extract the frequency
First step in the algorithm is centering the data by subtracting the sample mean. However, since we will compare with the true spectral density of the autoregression later, we will skip this step and compute the periodogram and spectral density of the true data.
The second step produces the initial estimate of the periodogram and
an estimated spectral density with a kernel \(K\) and a global bandwidth \(h\). The paper defines the periodogram to
be \[I(\omega) = \frac{1}{2\pi n} \Big|
\sum_{t=1}^n x_t e^{-it \omega} \Big|\] for \(\omega_j = 2\pi j/n, -n \leq j \leq n, -\pi \leq
\omega_j\leq \pi\).
The formula is, indeed, equivalent to our book’s definition of the
periodogram \[I(\omega_j) = |d(\omega_j)|^2 =
\frac{1}{n} \Big| \sum_{t=1}^n x_t e^{-2\pi i \omega_j t} \Big|\]
for \(j = 0,1,\dots, n-1\) since \(\omega = 2\pi j/n\), but in different
domains. Thus, we will the book’s code to compute the periodogram \(I\). Everything we does after this will be
in the frequency domain.
I = Mod(fft(x1)/sqrt(100))^2
I = I[1:floor(n/2)] # since I is symmetric, we will keep the first half
tsplot(freq, I, col=4)
Consequently, the paper chooses an inital bandwidth and a kernel to
calculate the estimated spectral density \[\hat{f}(\omega,h) = \frac{1}{nh} \sum_{k =-n}^n K
\Big(\frac{\omega-\omega_k}{h} \Big) I(\omega).\] where \(K(\cdot)\) is a given symmetric,
nonnegative function on the real values. It is noted that the inital
global bandwidth should not depend on \(\omega\). The paper chooses the Barlett
kernel with a global bandwidth of 0.05. In our exploration, we will work
with the Daniell kernel of spans
c(5,5). The next cell
gives a demonstration of the kernel weights and how we use the kernel to
get the smoothed spectral density.
(dm = kernel("modified.daniell", c(5,5))) # for a list
## mDaniell(5,5)
## coef[-10] = 0.0025
## coef[ -9] = 0.0100
## coef[ -8] = 0.0200
## coef[ -7] = 0.0300
## coef[ -6] = 0.0400
## coef[ -5] = 0.0500
## coef[ -4] = 0.0600
## coef[ -3] = 0.0700
## coef[ -2] = 0.0800
## coef[ -1] = 0.0900
## coef[ 0] = 0.0950
## coef[ 1] = 0.0900
## coef[ 2] = 0.0800
## coef[ 3] = 0.0700
## coef[ 4] = 0.0600
## coef[ 5] = 0.0500
## coef[ 6] = 0.0400
## coef[ 7] = 0.0300
## coef[ 8] = 0.0200
## coef[ 9] = 0.0100
## coef[ 10] = 0.0025
par(mfrow=1:2)
plot(kernel("daniell", 5), ylab=expression(h[~k]), col=4, panel.first=Grid(nxm=5))
plot(dm, ylab=expression(h[~k]), panel.first=Grid(), col=4) # for a plot
h = c(5,5) # initial bandwidth
kd = kernel("daniell", h)
f.hat = kernapply(I, kd, circular=TRUE)
tsplot(freq, f.hat, col=4)
Our next step is to estimate the residuals \[\hat{\epsilon} = \frac{I(\omega_k)}{\hat{f}(\omega_k,h_i)},\] for \(k = 1, \dots, n\). This results from the interpretation of the spectral estimator as the approximate multiplicative regression \[I(\omega_k) = \hat{f}(\omega_k) \cdot \epsilon_k.\] Then the paper recenters the data to avoid an additional bias in the resampling stage by rescaling the residuals \[\tilde{e}_k = \frac{\hat{\epsilon_k}}{\epsilon}, \text{ for } k=1, \dots, n, \text{ and } \epsilon = \frac{1}{n} \sum_{k=1}^n \hat{\epsilon}_k.\]
eps = I / f.hat
scaled_eps = eps / mean(eps)
mean(scaled_eps)
## [1] 1
The mean of the scaled residuals should be 1.
Now resample the scaled residuals and run boostrap estimate for a large number of times. The resampled periodogram is computed as \[I^*(\omega_k) = I^*(-\omega_k) = \hat{f}(w_k,g) \tilde{\epsilon_k}^*\] where \(g\) is a resampling bandwidth and \(\{\epsilon_k^*\}\) is a resample of the scaled residuals. We will pick \(g\) to be the same as the global bandwidth. Then the bootstrap kernel spectral density is calculated as \[\hat{f}(\omega,h,g) = \frac{1}{nh}\sum_{k=-n}^{n}K \Big(\frac{\omega-\omega_k}{h} \Big) I^*(\omega_k)\]
nboot = 1000
f.star.dist = c(c())
for (i in 1:nboot)
{
scaled_eps.star = sample(scaled_eps, replace=TRUE) # resample the scaled residuals
# new kernel
g = sample(1:5, 1)
kd.star = kernel("daniell", c(g,g))
# another estimated spectral density with the new kernel
f.hat = kernapply(I, kd.star, circular = TRUE)
# resampled periodogram
I.star = f.hat * scaled_eps.star
# bootstrap kernel spectral density
f.star = kernapply(I.star, kd, circular=TRUE)
f.star.dist[i] = list(f.star)
}
lwr = c()
upr = c()
for (i in 1:nspec)
{
f.star.wi = c()
for (j in 1:nboot)
{
f.star.wi[j] = f.star.dist[[j]][i]
}
# the 2.5 percentile spectral density at all frequencies
lwr[i] = quantile(f.star.wi, probs=c(.25, .75))[1]
# the 97.5 percentile spectral density at all frequencies
upr[i] = quantile(f.star.wi, probs=c(.25, .75))[2]
}
arma.spec(ar=c(0.5, -0.6, 0.3, -0.4, 0.2), main="Autoregression", ylim=c(0, max(upr,max(spec$spec))))
tsplot(freq, I, col=4)
lines(freq, upr, col=2)
lines(freq, lwr, col=2)
From the previous exploration, I build my first R package to bootstrap spectral density based on resampling the periodogram. You can find the complete package here: https://github.com/uyenle-gh/spec.boots
The following chunk installs the package into your environment.
Now we will apply our spec.boots function to produce a
confidence interval for the spectral density of the autoregression:
kd = kernel("daniell", c(5,5))
spec = spec.boots(x1, 1000, kd, demean=FALSE, alpha=0.95)
Most of the power is concentrated from frequency 0.13 to frequency 0.23, where the 95% confidence interval is above the median baseline. This allows us to conclude that the peak we observe around frequency 0.18 is different from the baseline and is statistically significant.
Another example where we know the exact frequencies is demonstrated below:
y1 = 2*cos(2*pi*1:100*6/100) + 3*sin(2*pi*1:100*6/100)
y2 = 4*cos(2*pi*1:100*10/100) + 5*sin(2*pi*1:100*10/100)
y3 = 6*cos(2*pi*1:100*40/100) + 7*sin(2*pi*1:100*40/100)
y = y1 + y2 + y3
y = ts(y)
spec.boots(y, 1000, kernel("daniell", c(2,2)))
## freq lwr upr
## 1 0.01 0.270223239 17.45207554
## 2 0.02 0.425149841 29.17483180
## 3 0.03 3.613353560 50.40524749
## 4 0.04 6.684144499 84.74913091
## 5 0.05 10.228294102 123.99912881
## 6 0.06 14.113791148 168.73174892
## 7 0.07 20.600388027 241.37377739
## 8 0.08 26.961871032 289.92339361
## 9 0.09 26.730406242 341.70273120
## 10 0.10 29.632080473 374.96964289
## 11 0.11 32.398298039 398.79006550
## 12 0.12 31.765245861 349.16376089
## 13 0.13 24.845023810 295.33638276
## 14 0.14 17.668298140 219.33958502
## 15 0.15 10.316856806 146.58849611
## 16 0.16 0.175864848 85.43525666
## 17 0.17 0.034660813 47.58141656
## 18 0.18 0.010259467 23.60354769
## 19 0.19 0.006406585 9.76800883
## 20 0.20 0.003795360 0.03046949
## 21 0.21 0.003587928 0.02883219
## 22 0.22 0.002863197 0.02676836
## 23 0.23 0.002688253 0.02466545
## 24 0.24 0.002506663 0.02304150
## 25 0.25 0.002309986 0.02119643
## 26 0.26 0.002289392 0.01988411
## 27 0.27 0.002019022 0.01818709
## 28 0.28 0.001860394 0.01632454
## 29 0.29 0.001854202 0.01600192
## 30 0.30 0.001750333 0.01493799
## 31 0.31 0.001617986 0.01457835
## 32 0.32 0.001594527 0.01409162
## 33 0.33 0.002667105 20.36772144
## 34 0.34 0.005072163 49.24809472
## 35 0.35 0.013794963 99.34851270
## 36 0.36 0.180623270 193.68520762
## 37 0.37 21.190621572 311.86491605
## 38 0.38 18.517373620 428.97680477
## 39 0.39 35.057150417 549.14836718
## 40 0.40 44.195429756 666.87758093
## 41 0.41 54.207926080 737.32322266
## 42 0.42 50.298037454 684.59572880
## 43 0.43 46.306856477 587.67380059
## 44 0.44 33.288276679 451.25448922
## 45 0.45 21.456315942 300.87576557
## 46 0.46 1.426250090 189.53943603
## 47 0.47 0.194874613 102.91956909
## 48 0.48 0.167327934 51.76333907
## 49 0.49 0.152880313 21.09192861
## 50 0.50 0.128516663 8.53530200
The spectral density plot of this example shows two peaks around frequency 0.1 and 0.4. The confidence interval corresponding to frequency 0.1 is slightly above the median baseline, so we suspect that this frequency only has a small effect.
In the next few chunks, we demonstrate how to use
spec.boots with real data and compare the result with those
produced by the mvspec command.
data("AirPassengers")
tsplot(AirPassengers, col=4, lwd=2)
mvspec(AirPassengers, col=4, lwd=2)
mvspec(AirPassengers, log='y', col=4, lwd=2)
mvspec(AirPassengers, log='y', spans=c(4,4), col=4, lwd=2)
spec.boots(AirPassengers, 1000, c(5,5))
## freq lwr upr
## 1 0.08333333 355.72568 4142.2920
## 2 0.16666667 420.54207 4576.5138
## 3 0.25000000 446.45384 4731.9612
## 4 0.33333333 422.74234 4440.5706
## 5 0.41666667 374.69175 3664.4454
## 6 0.50000000 370.13121 3068.4375
## 7 0.58333333 439.38391 3310.1799
## 8 0.66666667 603.41517 6251.2021
## 9 0.75000000 1010.98026 11819.4380
## 10 0.83333333 1672.14100 18909.4257
## 11 0.91666667 2470.73436 27376.6577
## 12 1.00000000 3079.55389 35076.7355
## 13 1.08333333 3237.30457 38125.7041
## 14 1.16666667 3126.18607 34922.3982
## 15 1.25000000 2541.43597 28261.3938
## 16 1.33333333 1819.00055 20618.4469
## 17 1.41666667 1149.11314 12820.0530
## 18 1.50000000 680.91683 7033.2147
## 19 1.58333333 470.55867 3831.4547
## 20 1.66666667 399.48993 3018.9933
## 21 1.75000000 473.26118 4337.7881
## 22 1.83333333 630.24126 6348.0246
## 23 1.91666667 851.76192 8489.9060
## 24 2.00000000 1096.04818 10263.8790
## 25 2.08333333 1064.68364 11140.4275
## 26 2.16666667 948.72250 10389.6724
## 27 2.25000000 676.99831 8208.0631
## 28 2.33333333 452.50726 5646.9596
## 29 2.41666667 261.12091 3266.7342
## 30 2.50000000 153.35273 1676.9989
## 31 2.58333333 89.47945 794.7784
## 32 2.66666667 69.14466 491.7707
## 33 2.75000000 69.27338 691.4979
## 34 2.83333333 95.04063 1011.3973
## 35 2.91666667 119.18149 1306.3992
## 36 3.00000000 140.21235 1499.2165
## 37 3.08333333 146.98432 1523.9207
## 38 3.16666667 133.40988 1400.2227
## 39 3.25000000 100.98750 1084.9823
## 40 3.33333333 70.90411 786.7359
## 41 3.41666667 43.53271 476.4370
## 42 3.50000000 31.06174 271.6371
## 43 3.58333333 28.12104 197.0744
## 44 3.66666667 33.69321 279.0043
## 45 3.75000000 48.24582 456.2280
## 46 3.83333333 67.65395 678.4795
## 47 3.91666667 90.14435 877.9652
## 48 4.00000000 111.81482 1113.2969
## 49 4.08333333 115.11076 1148.1218
## 50 4.16666667 104.20064 1066.9908
## 51 4.25000000 84.77034 877.5241
## 52 4.33333333 56.60093 665.5390
## 53 4.41666667 38.67698 434.6461
## 54 4.50000000 26.24454 254.1069
## 55 4.58333333 21.21035 157.6562
## 56 4.66666667 21.15284 174.4998
## 57 4.75000000 26.50933 271.4079
## 58 4.83333333 37.66147 424.9991
## 59 4.91666667 49.12580 557.1178
## 60 5.00000000 59.07278 665.9334
## 61 5.08333333 66.56840 710.2042
## 62 5.16666667 67.42700 700.9098
## 63 5.25000000 57.33906 621.9711
## 64 5.33333333 43.68773 484.2052
## 65 5.41666667 31.25137 300.4933
## 66 5.50000000 23.20855 186.3925
## 67 5.58333333 19.79099 133.5227
## 68 5.66666667 26.31766 238.2488
## 69 5.75000000 44.71610 579.7101
## 70 5.83333333 79.27307 1225.4516
## 71 5.91666667 148.83136 2084.3958
## 72 6.00000000 264.50198 3187.1431
For the AirPassengers dataset, it seems like both the
mvspec and spec.boots functions agree that
frequencies 1 and 2 are dominant. However, it is way more obvious
looking at the plot of spec.boots.
spec.boots(AirPassengers, 10, c(5,5))
## freq lwr upr
## 1 0.08333333 731.99173 3604.96187
## 2 0.16666667 716.04216 4498.92201
## 3 0.25000000 605.17727 4892.92272
## 4 0.33333333 484.18543 4286.74254
## 5 0.41666667 408.77396 3082.92833
## 6 0.50000000 469.72145 1934.98838
## 7 0.58333333 877.62054 2027.68736
## 8 0.66666667 964.02199 4084.04307
## 9 0.75000000 912.83499 8114.70850
## 10 0.83333333 1261.45085 13917.61787
## 11 0.91666667 2203.44924 21098.24015
## 12 1.00000000 2469.41840 28854.90366
## 13 1.08333333 2544.05172 34166.85402
## 14 1.16666667 2370.13739 34147.82709
## 15 1.25000000 2127.10441 31817.50394
## 16 1.33333333 1837.14736 25230.14253
## 17 1.41666667 1323.87855 16370.96073
## 18 1.50000000 882.92565 8788.68603
## 19 1.58333333 669.72833 3652.56033
## 20 1.66666667 425.39190 2304.67114
## 21 1.75000000 441.26783 2696.77237
## 22 1.83333333 943.47674 3284.44173
## 23 1.91666667 1499.53652 3652.39333
## 24 2.00000000 1884.16592 5099.91497
## 25 2.08333333 1875.04836 7413.34820
## 26 2.16666667 1535.66549 8914.19630
## 27 2.25000000 1243.78651 8557.34993
## 28 2.33333333 969.76814 6746.46504
## 29 2.41666667 471.10720 4321.30952
## 30 2.50000000 180.34888 2104.81837
## 31 2.58333333 91.92744 885.67222
## 32 2.66666667 86.47803 396.95258
## 33 2.75000000 117.95757 507.82314
## 34 2.83333333 131.60022 844.51740
## 35 2.91666667 142.94163 1234.72651
## 36 3.00000000 150.46923 1490.42414
## 37 3.08333333 162.82430 1473.97183
## 38 3.16666667 171.33678 1457.58254
## 39 3.25000000 170.63378 1221.12741
## 40 3.33333333 154.55544 857.35283
## 41 3.41666667 118.82683 497.88837
## 42 3.50000000 76.97959 248.99985
## 43 3.58333333 63.96134 169.27911
## 44 3.66666667 69.04767 224.12789
## 45 3.75000000 84.99141 397.98256
## 46 3.83333333 101.42293 652.75013
## 47 3.91666667 123.65424 850.83219
## 48 4.00000000 128.99579 878.63659
## 49 4.08333333 124.96755 867.35595
## 50 4.16666667 104.74997 785.47658
## 51 4.25000000 81.30024 617.65658
## 52 4.33333333 61.21470 473.33868
## 53 4.41666667 49.72262 318.54240
## 54 4.50000000 36.76382 195.37682
## 55 4.58333333 24.37189 111.65968
## 56 4.66666667 20.75320 162.23459
## 57 4.75000000 33.68745 261.83052
## 58 4.83333333 68.22190 360.08891
## 59 4.91666667 85.45678 405.19215
## 60 5.00000000 101.62507 422.70779
## 61 5.08333333 118.21557 427.27903
## 62 5.16666667 129.23131 362.59593
## 63 5.25000000 90.30335 270.69019
## 64 5.33333333 53.21480 176.58890
## 65 5.41666667 31.30471 122.13844
## 66 5.50000000 24.02900 88.96864
## 67 5.58333333 17.82785 70.17812
## 68 5.66666667 22.16462 166.49571
## 69 5.75000000 44.30087 345.86066
## 70 5.83333333 102.38449 760.90424
## 71 5.91666667 208.91635 1486.62636
## 72 6.00000000 410.45176 2526.63088
spec.boots(AirPassengers, 100, c(5,5))
## freq lwr upr
## 1 0.08333333 442.42580 3658.1634
## 2 0.16666667 481.83389 4457.9569
## 3 0.25000000 464.18697 4632.8774
## 4 0.33333333 382.23790 4036.2687
## 5 0.41666667 309.56109 3317.3710
## 6 0.50000000 415.41564 2960.4937
## 7 0.58333333 577.82897 3659.9509
## 8 0.66666667 659.41961 6687.0970
## 9 0.75000000 976.60470 11436.5491
## 10 0.83333333 1760.30341 19785.8486
## 11 0.91666667 2922.42179 28293.2361
## 12 1.00000000 4055.95350 35407.4041
## 13 1.08333333 4421.57554 35183.8524
## 14 1.16666667 3952.60320 32142.3481
## 15 1.25000000 3319.81450 27570.4327
## 16 1.33333333 2247.88169 20791.0557
## 17 1.41666667 1397.63082 14278.4184
## 18 1.50000000 746.19673 6998.9652
## 19 1.58333333 521.34407 3931.8703
## 20 1.66666667 457.15642 3356.4077
## 21 1.75000000 457.41300 4913.5418
## 22 1.83333333 606.77169 6773.6304
## 23 1.91666667 796.73626 8455.5933
## 24 2.00000000 930.02414 9482.5443
## 25 2.08333333 1050.95168 9635.2064
## 26 2.16666667 1079.67919 9337.3181
## 27 2.25000000 959.38954 8728.3364
## 28 2.33333333 664.49718 6536.2363
## 29 2.41666667 377.54695 3959.9556
## 30 2.50000000 168.19248 1712.9931
## 31 2.58333333 91.07768 746.9064
## 32 2.66666667 72.61043 488.8647
## 33 2.75000000 69.68855 720.5190
## 34 2.83333333 117.24534 1011.7487
## 35 2.91666667 142.20767 1317.6532
## 36 3.00000000 149.28078 1506.2503
## 37 3.08333333 147.26108 1490.7790
## 38 3.16666667 139.66640 1316.5922
## 39 3.25000000 115.01539 1063.6930
## 40 3.33333333 84.55511 757.6907
## 41 3.41666667 52.28924 431.5144
## 42 3.50000000 35.72845 254.5100
## 43 3.58333333 27.48639 183.2223
## 44 3.66666667 31.65797 268.9924
## 45 3.75000000 47.85126 426.0177
## 46 3.83333333 65.32690 632.5932
## 47 3.91666667 78.40531 825.6473
## 48 4.00000000 100.12381 941.9593
## 49 4.08333333 104.56446 1033.4007
## 50 4.16666667 93.32724 988.4380
## 51 4.25000000 80.13399 881.8727
## 52 4.33333333 63.21938 648.8008
## 53 4.41666667 46.48766 418.3601
## 54 4.50000000 29.58072 239.0013
## 55 4.58333333 27.05090 157.6391
## 56 4.66666667 24.53863 159.1314
## 57 4.75000000 32.84603 235.1024
## 58 4.83333333 43.57268 327.5970
## 59 4.91666667 53.95152 496.6233
## 60 5.00000000 61.92248 618.0584
## 61 5.08333333 58.70045 656.3071
## 62 5.16666667 57.21802 602.5993
## 63 5.25000000 44.78648 497.7824
## 64 5.33333333 31.71541 405.4326
## 65 5.41666667 27.21498 281.8897
## 66 5.50000000 20.82399 175.9036
## 67 5.58333333 17.42232 128.1382
## 68 5.66666667 22.21808 220.7756
## 69 5.75000000 39.19320 565.2491
## 70 5.83333333 102.80032 1002.6668
## 71 5.91666667 175.06464 1584.4272
## 72 6.00000000 306.94794 2488.0804
spec.boots(AirPassengers, 1000, 2)
## freq lwr upr
## 1 0.08333333 341.444302 4277.54939
## 2 0.16666667 740.944938 6529.71298
## 3 0.25000000 806.370292 6064.23136
## 4 0.33333333 652.884823 4420.27432
## 5 0.41666667 340.877207 2531.04294
## 6 0.50000000 144.015147 1031.45656
## 7 0.58333333 75.744429 501.75349
## 8 0.66666667 73.678918 500.72492
## 9 0.75000000 79.457507 662.10614
## 10 0.83333333 207.642851 1812.37740
## 11 0.91666667 1082.558324 15001.96035
## 12 1.00000000 4922.714994 47052.32894
## 13 1.08333333 8839.739796 63608.47780
## 14 1.16666667 5591.124442 48459.96770
## 15 1.25000000 1360.666402 16253.96397
## 16 1.33333333 381.701885 2931.43967
## 17 1.41666667 226.720827 1404.56653
## 18 1.50000000 221.561571 1403.54418
## 19 1.58333333 213.597399 1481.59931
## 20 1.66666667 177.011833 1164.71857
## 21 1.75000000 121.873145 832.56759
## 22 1.83333333 145.118377 1008.82786
## 23 1.91666667 455.061439 4932.88245
## 24 2.00000000 1521.120382 14380.27723
## 25 2.08333333 2697.238630 19114.68565
## 26 2.16666667 1191.207706 13605.74107
## 27 2.25000000 244.195974 4039.73821
## 28 2.33333333 43.900488 370.92427
## 29 2.41666667 16.472658 117.60308
## 30 2.50000000 9.120546 59.26724
## 31 2.58333333 11.121379 96.11007
## 32 2.66666667 26.452875 196.90107
## 33 2.75000000 38.917312 263.12387
## 34 2.83333333 52.924612 325.61515
## 35 2.91666667 106.645849 846.37170
## 36 3.00000000 237.579059 2057.96118
## 37 3.08333333 305.460200 2529.64834
## 38 3.16666667 163.555070 1703.30822
## 39 3.25000000 34.845390 517.60535
## 40 3.33333333 8.607293 60.39535
## 41 3.41666667 6.645641 41.13489
## 42 3.50000000 7.519997 52.66790
## 43 3.58333333 12.854080 88.38006
## 44 3.66666667 17.553135 119.04435
## 45 3.75000000 18.744935 121.74053
## 46 3.83333333 24.158900 164.85127
## 47 3.91666667 64.058979 593.15868
## 48 4.00000000 167.670579 1472.08313
## 49 4.08333333 225.920178 1809.24444
## 50 4.16666667 161.465362 1280.96584
## 51 4.25000000 67.777985 542.57270
## 52 4.33333333 26.190474 221.07157
## 53 4.41666667 11.019836 87.93691
## 54 4.50000000 7.034112 44.25267
## 55 4.58333333 8.442239 63.71849
## 56 4.66666667 8.418808 70.94649
## 57 4.75000000 8.162076 52.16154
## 58 4.83333333 7.880947 64.42987
## 59 4.91666667 30.507102 297.38006
## 60 5.00000000 94.033291 858.56049
## 61 5.08333333 151.174808 1132.74902
## 62 5.16666667 124.259334 884.01798
## 63 5.25000000 48.292974 383.43208
## 64 5.33333333 15.770809 139.53200
## 65 5.41666667 8.875523 63.97610
## 66 5.50000000 7.971735 51.30126
## 67 5.58333333 8.868466 64.49239
## 68 5.66666667 9.305446 62.72120
## 69 5.75000000 7.541146 51.22882
## 70 5.83333333 6.760781 47.46511
## 71 5.91666667 5.802388 37.91271
## 72 6.00000000 31.210300 1081.91102
spec.boots(AirPassengers, 1000, 8)
## freq lwr upr
## 1 0.08333333 245.72442 3325.7598
## 2 0.16666667 284.96873 3684.2401
## 3 0.25000000 299.99114 3704.7681
## 4 0.33333333 316.17489 3582.7041
## 5 0.41666667 335.45840 4659.5291
## 6 0.50000000 410.30992 7920.3407
## 7 0.58333333 567.40403 12250.2817
## 8 0.66666667 837.07457 16240.4331
## 9 0.75000000 1080.28118 19260.9771
## 10 0.83333333 1358.12008 21988.4764
## 11 0.91666667 1646.66765 23984.4010
## 12 1.00000000 1990.82808 25968.8320
## 13 1.08333333 2258.46920 27294.4957
## 14 1.16666667 2064.97024 25949.3912
## 15 1.25000000 1621.02419 23858.2867
## 16 1.33333333 1234.87564 21728.3801
## 17 1.41666667 1021.76312 19600.5633
## 18 1.50000000 885.26480 17046.5163
## 19 1.58333333 705.89595 13988.9355
## 20 1.66666667 614.27214 9512.9975
## 21 1.75000000 537.29905 7190.3319
## 22 1.83333333 537.91562 7273.0448
## 23 1.91666667 595.96259 7896.4912
## 24 2.00000000 693.29845 8293.8942
## 25 2.08333333 698.88043 8448.2797
## 26 2.16666667 609.63972 8027.4040
## 27 2.25000000 458.91375 7268.5592
## 28 2.33333333 353.80351 6635.3483
## 29 2.41666667 277.37737 5751.2527
## 30 2.50000000 198.08170 4636.9654
## 31 2.58333333 149.80058 3453.2663
## 32 2.66666667 103.98952 2156.2062
## 33 2.75000000 82.94587 1199.6533
## 34 2.83333333 75.86704 1024.3543
## 35 2.91666667 84.88908 1090.1318
## 36 3.00000000 94.79212 1184.2171
## 37 3.08333333 95.56500 1199.9821
## 38 3.16666667 88.34112 1159.7768
## 39 3.25000000 67.24408 1032.6274
## 40 3.33333333 56.39605 935.4418
## 41 3.41666667 49.02389 795.5725
## 42 3.50000000 49.41510 704.6571
## 43 3.58333333 43.23553 634.0099
## 44 3.66666667 43.22042 601.4500
## 45 3.75000000 42.68042 653.4867
## 46 3.83333333 44.47090 712.7974
## 47 3.91666667 56.59589 759.1904
## 48 4.00000000 69.04058 827.8668
## 49 4.08333333 74.32640 872.2511
## 50 4.16666667 61.02890 826.9887
## 51 4.25000000 53.83024 750.9795
## 52 4.33333333 45.38397 697.4350
## 53 4.41666667 38.20250 601.1231
## 54 4.50000000 34.62182 536.8590
## 55 4.58333333 32.89864 487.6930
## 56 4.66666667 29.34218 399.1542
## 57 4.75000000 27.17421 400.1414
## 58 4.83333333 32.64471 437.2461
## 59 4.91666667 37.98641 477.6107
## 60 5.00000000 42.13411 513.5739
## 61 5.08333333 48.06889 549.5802
## 62 5.16666667 42.95257 540.1455
## 63 5.25000000 35.20777 505.6799
## 64 5.33333333 29.73180 445.8542
## 65 5.41666667 23.19520 384.4005
## 66 5.50000000 25.34589 393.4845
## 67 5.58333333 29.07637 621.9836
## 68 5.66666667 36.81357 1137.1677
## 69 5.75000000 65.48339 1568.7786
## 70 5.83333333 99.99604 2034.1136
## 71 5.91666667 131.41770 2536.8195
## 72 6.00000000 169.15290 2965.7731
spec.boots(AirPassengers, 1000, c(2,2))
## freq lwr upr
## 1 0.08333333 443.549063 4239.44806
## 2 0.16666667 639.606722 6087.67741
## 3 0.25000000 644.087326 6098.44500
## 4 0.33333333 554.398750 4904.26068
## 5 0.41666667 373.617993 3031.48960
## 6 0.50000000 219.755387 1588.61827
## 7 0.58333333 123.951552 788.07257
## 8 0.66666667 105.811598 663.83454
## 9 0.75000000 229.711747 1977.77862
## 10 0.83333333 737.620155 8070.26011
## 11 0.91666667 2161.952943 24665.15946
## 12 1.00000000 4301.453835 45344.60117
## 13 1.08333333 5548.850714 54523.10612
## 14 1.16666667 4203.123382 44265.22524
## 15 1.25000000 2384.658099 24909.19166
## 16 1.33333333 872.442815 9244.49444
## 17 1.41666667 390.453377 2707.94767
## 18 1.50000000 246.347147 1603.82899
## 19 1.58333333 193.706854 1444.01550
## 20 1.66666667 178.262057 1262.90474
## 21 1.75000000 218.399711 1256.50583
## 22 1.83333333 401.149342 2931.11612
## 23 1.91666667 796.698669 7419.03105
## 24 2.00000000 1412.073153 13691.30298
## 25 2.08333333 1576.930956 16664.18008
## 26 2.16666667 1309.456827 13360.59793
## 27 2.25000000 653.242295 6906.60356
## 28 2.33333333 203.649287 2159.11667
## 29 2.41666667 44.476145 474.62204
## 30 2.50000000 19.275123 110.81347
## 31 2.58333333 21.269207 137.12289
## 32 2.66666667 28.043769 213.92030
## 33 2.75000000 45.508084 319.42670
## 34 2.83333333 81.196553 562.45035
## 35 2.91666667 151.853571 1163.12013
## 36 3.00000000 212.251386 1955.79692
## 37 3.08333333 246.104074 2284.08886
## 38 3.16666667 198.943810 1717.12585
## 39 3.25000000 88.157813 916.19108
## 40 3.33333333 29.009348 312.28033
## 41 3.41666667 12.573255 85.70476
## 42 3.50000000 9.745460 64.81780
## 43 3.58333333 11.992911 89.22978
## 44 3.66666667 15.662776 115.46497
## 45 3.75000000 24.858639 163.79220
## 46 3.83333333 46.407357 373.01514
## 47 3.91666667 92.111412 817.71123
## 48 4.00000000 148.079661 1333.61939
## 49 4.08333333 164.037568 1581.37678
## 50 4.16666667 140.535882 1238.32119
## 51 4.25000000 86.878977 715.37408
## 52 4.33333333 39.267808 321.10179
## 53 4.41666667 19.291789 138.96250
## 54 4.50000000 11.166202 73.87106
## 55 4.58333333 9.810999 65.09039
## 56 4.66666667 9.159427 66.09173
## 57 4.75000000 11.721640 75.79076
## 58 4.83333333 22.563940 182.73917
## 59 4.91666667 45.961964 465.10797
## 60 5.00000000 86.595899 761.39734
## 61 5.08333333 112.591694 975.62279
## 62 5.16666667 97.259228 841.37928
## 63 5.25000000 62.207332 532.61801
## 64 5.33333333 32.544075 256.46344
## 65 5.41666667 16.951300 109.01330
## 66 5.50000000 11.179473 67.02855
## 67 5.58333333 8.869169 64.90745
## 68 5.66666667 8.310764 63.58237
## 69 5.75000000 7.833013 57.44793
## 70 5.83333333 14.717061 123.50386
## 71 5.91666667 51.737456 604.75272
## 72 6.00000000 192.261129 2005.72210
spec.boots(AirPassengers, 1000, c(8,8))
## freq lwr upr
## 1 0.08333333 260.39016 3335.9897
## 2 0.16666667 300.09014 3976.1469
## 3 0.25000000 376.92404 4780.9006
## 4 0.33333333 487.57725 5961.1997
## 5 0.41666667 596.99555 7736.3405
## 6 0.50000000 750.96404 9610.3397
## 7 0.58333333 930.72594 11503.3144
## 8 0.66666667 1047.39265 14103.5523
## 9 0.75000000 1164.85491 16454.2744
## 10 0.83333333 1295.45320 18549.0080
## 11 0.91666667 1398.17112 19763.5414
## 12 1.00000000 1477.67391 21015.2692
## 13 1.08333333 1518.95767 21303.2890
## 14 1.16666667 1481.74572 21302.5647
## 15 1.25000000 1438.46224 19926.3435
## 16 1.33333333 1373.35037 18468.9810
## 17 1.41666667 1398.53224 16859.1364
## 18 1.50000000 1229.68647 14905.7020
## 19 1.58333333 1077.74885 12563.6650
## 20 1.66666667 965.14142 11013.6439
## 21 1.75000000 864.62525 9578.9880
## 22 1.83333333 768.65372 8639.2735
## 23 1.91666667 679.24258 7790.0533
## 24 2.00000000 623.58718 7181.6051
## 25 2.08333333 553.80602 6910.1308
## 26 2.16666667 493.66580 6615.9237
## 27 2.25000000 438.28011 6044.7737
## 28 2.33333333 389.05269 5306.5804
## 29 2.41666667 337.77389 4586.1084
## 30 2.50000000 299.30630 3994.4885
## 31 2.58333333 263.90647 3316.3387
## 32 2.66666667 220.28405 2777.8399
## 33 2.75000000 190.37876 2191.1656
## 34 2.83333333 157.28232 1790.9803
## 35 2.91666667 134.12916 1447.9269
## 36 3.00000000 110.46740 1293.5977
## 37 3.08333333 94.35345 1139.7722
## 38 3.16666667 83.97886 1043.8100
## 39 3.25000000 78.70124 943.8301
## 40 3.33333333 68.24512 853.2943
## 41 3.41666667 65.04799 785.2474
## 42 3.50000000 63.36960 719.6929
## 43 3.58333333 60.81841 686.4710
## 44 3.66666667 59.65766 690.6790
## 45 3.75000000 58.47671 669.5236
## 46 3.83333333 57.04074 690.4504
## 47 3.91666667 56.76413 673.7383
## 48 4.00000000 55.00133 667.7482
## 49 4.08333333 55.47400 694.3172
## 50 4.16666667 54.57006 669.9030
## 51 4.25000000 51.32294 669.5381
## 52 4.33333333 50.10691 610.1789
## 53 4.41666667 49.33619 596.9437
## 54 4.50000000 47.08598 539.1673
## 55 4.58333333 46.90944 481.2551
## 56 4.66666667 42.48079 444.6242
## 57 4.75000000 41.26729 436.1385
## 58 4.83333333 39.11246 429.7523
## 59 4.91666667 36.78001 427.8580
## 60 5.00000000 37.68327 435.7236
## 61 5.08333333 37.90292 454.5900
## 62 5.16666667 41.28003 459.8263
## 63 5.25000000 44.81580 473.7121
## 64 5.33333333 50.50933 534.8685
## 65 5.41666667 60.86444 653.4778
## 66 5.50000000 69.67280 844.1693
## 67 5.58333333 84.50423 1064.5311
## 68 5.66666667 102.98324 1364.0618
## 69 5.75000000 120.44491 1630.9000
## 70 5.83333333 145.19924 1935.1037
## 71 5.91666667 176.83469 2424.5463
## 72 6.00000000 208.97888 2854.7019
One thing to note here is that the result heavily depends on the global bandwidth.
data(Inflation)
tsplot(Inflation$CPIPctDiff, col=4, lwd=2)
mvspec(Inflation$CPIPctDiff, col=4, lwd=2)
mvspec(Inflation$CPIPctDiff, col=4, lwd=2, log="y")
mvspec(Inflation$CPIPctDiff, col=4, lwd=2, log="y", spans=c(4,4))
kd = kernel("daniell", c(5,5))
spec.boots(ts(Inflation$CPIPctDiff) , 1000, kd)
## freq lwr upr
## 1 0.01041667 0.04748392 0.16181563
## 2 0.02083333 0.05325450 0.17933588
## 3 0.03125000 0.05952752 0.19781128
## 4 0.04166667 0.06512552 0.21691661
## 5 0.05208333 0.07069815 0.23131778
## 6 0.06250000 0.07539897 0.24772430
## 7 0.07291667 0.08092185 0.26408992
## 8 0.08333333 0.08412068 0.27315913
## 9 0.09375000 0.08755368 0.28310885
## 10 0.10416667 0.08896712 0.29471411
## 11 0.11458333 0.09085065 0.29900890
## 12 0.12500000 0.08983223 0.30410073
## 13 0.13541667 0.09008047 0.30049803
## 14 0.14583333 0.08916984 0.29745182
## 15 0.15625000 0.08836796 0.28632948
## 16 0.16666667 0.08530142 0.27067392
## 17 0.17708333 0.08273448 0.25578083
## 18 0.18750000 0.07723313 0.24513634
## 19 0.19791667 0.07279217 0.22922028
## 20 0.20833333 0.06810377 0.20846959
## 21 0.21875000 0.06269869 0.19259222
## 22 0.22916667 0.05700180 0.17470333
## 23 0.23958333 0.05188547 0.16040058
## 24 0.25000000 0.04666218 0.14429562
## 25 0.26041667 0.04176271 0.12975493
## 26 0.27083333 0.03780405 0.11510366
## 27 0.28125000 0.03403855 0.10417829
## 28 0.29166667 0.03049161 0.09396414
## 29 0.30208333 0.02748123 0.08420837
## 30 0.31250000 0.02493162 0.07543659
## 31 0.32291667 0.02344686 0.06760640
## 32 0.33333333 0.02176075 0.06222577
## 33 0.34375000 0.02068953 0.05698164
## 34 0.35416667 0.01947760 0.05519533
## 35 0.36458333 0.01872989 0.05290752
## 36 0.37500000 0.01802600 0.05252007
## 37 0.38541667 0.01809865 0.05224218
## 38 0.39583333 0.01827226 0.05351587
## 39 0.40625000 0.01875471 0.05603375
## 40 0.41666667 0.01991227 0.05823241
## 41 0.42708333 0.02156429 0.06235975
## 42 0.43750000 0.02362367 0.06838759
## 43 0.44791667 0.02541298 0.07669854
## 44 0.45833333 0.02821586 0.08744856
## 45 0.46875000 0.03065814 0.09791672
## 46 0.47916667 0.03353582 0.11233245
## 47 0.48958333 0.03759505 0.12861959
## 48 0.50000000 0.04204738 0.14335926
The spectral density of the CPI returns peaks around frequency 0.1
for both datasets. Here, it is easier to tell which frequencies matter
since we have a baseline to rely on. One drawback to this approach, as
you can also tell from the two plots of spec.boots and
mvspec, is that some of the peaks are flattened and spread
out to others due to averaging.